Cosmological dynamics of fourth order gravity: A compact view 
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We construct a compact phase space for flat FLRW spacetimes with standard matter described 
by a perfect fluid with a barotropic equation of state for general f(R) theories of gravity, subject 
to certain conditions on the function /. We then use this framework to study the behaviour of the 
phase space of Universes with a non-negative Ricci scalar in R + aR n gravity. We find a number 
of interesting cosmological evolutions which include the possibility of an initial unstable power-law 
inflationary point, followed by a curvature fluid dominated phase mimicking standard radiation, 
then passing through a standard matter (CDM) era and ultimately evolving asymptotically towards 
a de-Sitter-like late-time accelerated phase. 
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INTRODUCTION 

The ACDM (or Concordance) Model p} is one of the 
greatest successes of General Relativity. It reproduces 
beautifully all the main observational results e.g., the 
dimming of type la Supernovae Q, Cosmic Microwave 
Background Radiation (CMBR) anisotropies Q, Large 
Scale Structure formation 0] , baryon oscillations [B[ and 
weak lensing Q). Unfortunately, this model is also af- 
fected by significant fine-tuning problems related to the 
vacuum energy scale and this has led to a considerable 
amount of effort throughly exploring other viable theo- 
retical schemes. 

Currently, one of the most popular alternatives to the 
ACDM model is based on gravitational actions which are 
non-linear in the Ricci curvature R and/or contain terms 
involving combinations of derivatives of R: the so called 
f(R) theories of gravity [9l4l2|. Such models first be- 
came popular in the 1980's because it was shown that 
they can be derived from fundamental physical theories 
(Like M-theory) and naturally admit a phase of acceler- 
ated expansion, which could be associated with an early 
universe inflationary phase Q. The fact that the phe- 
nomenology of Dark Energy requires the presence of a 
similar phase (although only a late time - low energy 
one) has recently revived interest in these models. In 
particular, the idea that Dark Energy may have a geo- 
metrical origin, i.e., that there is a connection between 
Dark Energy and a non-standard behavior of gravitation 
on cosmological scales is now a very active area of re- 
search (see for example 13- 19j|). 

Unfortunately efforts to obtain an understanding of the 
physics of these theories are hampered by the complexity 
of the fourth-order field equations, making it difficult to 
obtain both exact and numerical solutions, which can be 
compared with observations. Recently however, progress 



has been made in resolving these issues using a number of 
useful techniques. One such method, based on the theory 
of dynamical systems [5(| has proven to be very successful 
in providing a simple way of obtaining exact solutions 
and a (qualitative) description of the global dynamics of 
these models (22|. The dynamical systems analysis has 
up to now approached this problem in the conventional 
way, by first exploring the finite equilibrium points and 
then computing the asymptotic behavior of the phase- 
space using the classical method of Poincare projections. 

In this paper we develop an alternative scheme which 
involves compactifying the phase space for general f(R) 
theories of gravity, subject to certain conditions on the 
function /. We then use this framework to study the 
behaviour of the phase space of Universes with a non- 
negative Ricci scalar in R + aR n gravity. We find a 
number of interesting cosmological evolutions, which in- 
clude the possibility, at least in principle, where the 
Universe begins close to at an unstable power-law infla- 
tionary equilibrium point, then evolves towards a curva- 
ture fluid dominated phase where the effective equation 
of state mimics standard radiation with w *~ 1/3 (we 
will refer to such phases as radiation-like), then passes 
through a standard matter (CDM) era and ultimately 
evolves asymptotically towards a de-Sitter-like late-time 
accelerated phase. 

We also show that as n —> 0, all the fixed points 
that approach the ACDM subspace of the complete state 
space of R + aR n gravity, are unstable. This implies 
that the behavior of the solutions of a fourth order the- 
ory which is close to ACDM may be completely different 
from those of ACDM and one needs to do a careful anal- 
ysis of the solutions rather than a priori assuming any 
global behavior of the trajectories. 
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THE FIELD EQUATIONS FOR FOURTH ORDER 
GRAVITY 



The natural extension of standard General Relativity is 
to consider a Lagrangian that contains curvature invari- 
ants of higher than linear order. In fact, renormalization 
of quantum field theory suggests that adding such terms 
to the standard gravitational action appears to be neces- 
sary [24| to give a first approximation to some quantized 
theory of gravity. The quadratic Lagrangians are at the 
first level of such modifications and have been studied ex- 
tensively over true past two decades. The four possible 
second-order curvature invariants are 
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where e lklm is completely antisymmetric 4-volume ele- 
ment and R, R ab , Rabcd are the Ricci scalar, Ricci tensor 
and Riemann tensor respectively. However, for homoge- 
neous and isotropic spacetimes, because of the following 
identities 23] 
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it follows that the general fourth-order Lagrangian for a 
these highly symmetric spacetimes contain only powers 
of R and we can write the action as 



A =\j d i x^—g[f{R) + 2C m ] 



(5) 



where C m represents the matter contribution. 

Varying the action with respect to the metric gives the 
following field equations: 

f'G ab = r^ + i(/-i?/').g Qb + V h V a /'-.g ah V c V c /', (6) 

where /' denotes the derivative of the function / w.r.t. 
the Ricci scalar and T"J is the matter stress energy tensor 
defined by 

T2 = p m u a u b + p m h ab + q™u b + ql n u a + 7r™ b . (7) 

Here u a is the direction of a timelike observer, h ab is 
the projected metric on the 3-space perpendicular to u a . 
Also fi m , p m , q m and ir" l b denotes the standard matter 
density, pressure, heat flux and anisotropic stress respec- 
tively. Equations ^1 reduce to the standard Einstein 
field equations when f(R) = R. 

For the homogeneous and isotropic spacetimes with 
vanishing 3-curvature and barotropic perfect fluid as the 
standard matter source with equation of state p = ujp, 
the independent field equations for general f(R) gravity 
are as follows. 



• The Raychaudhuri equation: 
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where O is the volume expansion of the matter flow 
lines u a and p is the standard matter density. 



• The Friedmann equation: 
/' 2 2/ 

and 



p>2_3p 3 3/ (/') 



• Conservation of standard matter. 

p= -Q(l + u>)p. 



(10) 



Combining the Raychaudhuri and Friedman equations, 
we obtain: 



R = 29 + -e 2 . 



(11) 



COMPACT PHASE SPACE FOR POSITIVE 
RICCI SCALAR UNIVERSE: 

In this paper we will study the dynamics of Friedmann - 
Lemaitre- Robertson- Walker (FLRW) models only in 
the sector R > 0. This is because the sector R < 
is not of much physical interest and also, as we shall see 
later, the sectors R > and R < are connected by 
the invariant sub- manifold R = 0, making the physically 
interesting dynamics completely confined to the sector 
R > 0. Also, we consider the 3-curvature to be vanishing, 
which is an invariant sub-manifold by itself. As required 
by the no -ghost condition we also assume /' > 0. 

To compactify the phase space we rewrite the Fried- 
mann equation ([9]) in the following form: 



where 
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We can now define the following set of normalized vari- 
ables: 



sen 

2f'D 



3 / 



2 f'D 2 

_ 3 R 
Z ~2Tfi 



f'D 2 

Q = %- (14) 



3 



To guarantee that the propagation equations for these 
compact variables will result in a dimensionless dynami- 
cal system, we need to define a new time variable r, such 
that 
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For r to be a monotonously increasing time variable, a 
normalization D is chosen such that it is strictly positive 
at all times. It is clear by construction that when = 
the normalized dynamical variables as well as the time 
variable are well defined, thus this normalization allows 
the study of general static, re-collapsing and bouncing 
solutions. 

From the Friedmann equation we obtain the following 
constraints, 
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The first constraint comes directly from Friedmann equa- 
tion, while the second one arises from the definition of 
the normalization parameter D. According to these con- 
straints and considering R > 0, p > and /' > 0, we see 
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that the above dynamical variables have to be defined in 
the following ranges, 

< Q m < 1 , < z < 1 , -1 < x < 1 

-2 < Q < 2; < y < 1 , (17) 

making the complete phase space compact. Also since 
the variable Q is a normalized Hubble parameter, the cos- 
mological solutions will naturally include both expanding 
and collapsing as well as static solutions and these two 
sets of solutions are connected via the non-invariant sub- 
set Q = 0. 



THE PROPAGATION EQUATIONS 

An autonomous system, which is equivalent to cosmo- 
logical equations (jSlTTTj) can be derived by differentiating 
the compact variables (|14p . with respect to r and using 
(|51 ITU|) . The dimensionality of the resultant system can 
then be reduced by using the two constraints (IT6l) . By 
eliminating the dynamical variables Q m and y, we obtain 
the following 3-dimensional effective autonomous system: 



x' = - ( - 3(1 + u) - (1 + 3w)a; 4 - 4Q 2 (-1 + x 2 ) + (1 + 3u)z 
6 V 
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+ x(5 + 3cj - (1 + 3lj)x 2 - T - 2 - 3n(l + u) + 2z + 3n(l + u)(z + (Q + x) 2 ) ^ 
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+ x(5 + 3w - (1 + 3oj)x 2 - T - 2 - 3n(l + w) + 2z + 3n(l + w)((Q + a;) 2 + z) ^ , 



where T = f'/Rf". In general, the system is not closed 
unless r is expressed in terms of the dynamical variables 
(fT4"]) . For example, in the case of R + aR n , we have 



n(y — z) n[{q + x) 2 + z — 1] 

Thus, the above system defines the dynamics of all well 
defined f(R) theories for which f'/Rf" is invertible in 
terms of the dynamical variables. From equations (|18[) 
we can see that z = is an invariant sub-manifold, and 
in the z — 2-surface the line Q = is an invariant 
subset. Since z = corresponds to R — 0, we obtain an 
important result: 

For all well defined functions f(R), with f > and 
f'/Rf" invertible in terms of the dynamical variables de- 
fined by |J^[ ), a FLRW universes with non - negative Ricci 
Scalar, continues to be so, both in the future and in the 



past. Also an R — universe can never undergo a bounce 
in the future or the past. 



In the next section we will fix the function / to be the 
class of theories f(R) = R+aR n and study the dynamics 
of the flat FLRW universes and their stability for those 
theories. In order to study the stability of the fixed points 
of the dynamical systems (fT8"|) , we will use the very well 
known techniques, which involve linearizing the dynam- 
ical equations around the equilibrium points and then 
finding the eigenvalues of the linearization matrix (the 
Jacobian) at the equilibrium points. If the Jacobian is 
well defined, then they can be classified according to the 
sign of the real part of eigenvalues as attractors, repellers 
and saddle points. 
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THE FIXED POINTS AND EXACT SOLUTIONS 
FOR R + aR n GRAVITY 



As we have seen from the equation p^|) , f'/Rf" is 
invertible in terms of the dynamical variables for f(R) — 
R + aR n . It is interesting to note that the constant 
'n' couples to the dynamical equations (fl"8|) only via the 
quantity V and the constant a does not couple to the 
equations at all. Hence all the fixed point of the system 
are necessarily independent of a. 



The coordinates of the fixed points are shown in Table 
1. Note that each fixed point has an expanding (Q > 0) 
and a collapsing [Q < 0) version as indicated by the sub- 
scripts (+; — ) respectively. Also some points only occur 
in the compact state space defined by (IT71) for certain 
ranges of n. The occurrence of the fixed points outside 
the compact region for specific n and uj means that the 
constraints (fT6|) are not satisfied and consequently these 
fixed points are not physical for these values of n and ui). 
Fixed points that are not physical for these values of n 
and ui have been excluded from the analysis. 



Table 1. Coordinates of the equilibrium points for R + aR n -gravity. We will not explicitly state the expressions 
for s,gl,...,g4 and fl,...,f4, which are rational functions of n and u, however we give them at the following link [2lj ], 



Fixed points 



A ± 
B 

C 
D ± 

E± 

F± 
G± 

I± 

L± 

N ± 



Coordinates (x, il, z, Q) 



(1,0,0, ±2) 
(±1,0,0,0) 



(- 

I 

( 



V3 + 12c; + 9a; 2 

1 + 3cj 
- 3a; 4(3w - 



2) 



l + 3u; 
,0,± 



0,0) 



3(w-l)' 9(w-l) 2 ' ' 3(w-l)' 
(0,0,1,±-J=) 

(/i (n, u),gi (n, u),h (n, id), m (n, lu)) 
(/ 2 (n, ui),g 2 (n, ui),l 2 (n, uj), n 2 {n, u)) 



(Mn),g a (n),l 3 (n),n 3 (n)) 

(fA{n,u),g A (n,u),l A (n,u)),n A (n,u))) 
(0,^,±f ) 
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o ((n - 2)t - t ) 



-l + 3n-2n 2 
-2 + n 
2n 



a {Zt{l + u)~to)3( l +u) 
a (2t-to) 2/3 



By looking at the coordinates of the fixed points in Table 
1, we can distinguish two classes; the first corresponds to 
points with coordinates that are independent of n, which 
means that these points are common to all f(R) theories. 
This class contains the fixed points A±,B,C±,D±,E± 
and N± and they all lie on the boundary of the compact 
region except for the point N. 

In the non- compact analysis developed in [26| . non 
of these boundary points appear. Furthermore, even 
though N± is not a boundary point, it doesn't appear 
in [26|, because of it's special location in the phase space 
- it lies exactly on the intersection of the plane x = and 
the surface z — y = 1 — (Q + x) 2 . In this case one has 
to take the limit of T carefully as one approaches this 
point and the standard techniques of finding fixed points 
breaks down for this case. 

The other class contains fixed points with coordinates 
that depend on n and w; this class contains the three 
points L± , I± and F± . F± is the only boundary point 
and it lies in the invariant sub-manifold z = 0. The 
expanding versions of the points L± and I± correspond 
to the equally labeled finite points in 26(. The point H 



in [26[ enters the compact sector, which we consider in 
this paper only when n = (1 + \/3)/2 and for this value 
of n it merges with the point /. All the other points that 
appear in the above mentioned reference do not appear 
in the sector we are studying in this paper. 
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Table 2. The stability of the fixed points for ui = 0; 1/3. 



Fixed point 



Physical range 



T 



Stability 



Equation of state 



LJ = 



OJ = 



57, 



Vn 
Vn 
Vn 
Vn 
Vn 
Vn 

^(0,i + — J 
n G (-, 1) and n > 5/4 



.3 4 

" G( 4'7 
Vn 

Vn 



Vn 
Vn 
Vn 
Vn 
Vn 
Vn 



37, 



l vTf. 



ne (0,- 

1 8 ° 

n G (-, 1) and n > 5/4 



«e (l,\/2) 

Vn 
Vn 



Attractor 

Repeller 

Attractor 

Saddle 

Repeller 

Attractor for n € (0, 2) 
Saddle 

Saddle for n€ (1/2,1) 

Attractor for n = 5/4 
Saddle for n G (5/4,2) 
Attractor for n > 2 

Saddle for n G (1/2,1) 

Repeller for tj = 5/4 
Saddle for n G (5/4,2) 
Repeller for 71 > 2 

Saddle 

Spiral + 
Spiral - 



Attractor 

Repeller 

Attractor 

Saddle 

Repeller 

Attractor for n e (0, 2) 
Saddle 

Saddle for n G (1/2,1) 

Attractor for n = 5/4 
Saddle for n G (5/4,2) 
Attractor for n > 2 

Saddle for n G (1/2,1) 

Repeller for n = 5/4 
Saddle for n G (5/4,2) 
Repeller for n > 2 

Saddle 

Spiral + 
Spiral - 



Exact solutions 

The exact solutions at the fixed points are also sum- 
marized in Table 1 and the stability analysis for the dust 
and radiation cases are summarized in Table 2. First 
we discuss the static solutions. From definition (fl8|. 
Q = => 6 = 0, so any fixed point that lies on the 
surface Q = represents a static Universe. By looking 
at the coordinates of the fixed points in Table 1, we can 
see that the point B is static for all values of n and to. 
The point I± is static only for n = 1/2 and n = 1 and 
we find that for these values of n the point I± represent 
an unstable saddle point. 

We now proceed to find the exact solutions for the scale 
factor at the non- static fixed points. The expansion rate 

da 

and the deceleration parameter q — are related by 

a 

the Raychaudhuri equation, 

e = -i(i + g )e 2 . (20) 

If we know the value of the deceleration parameter qi at 
some fixed point i, we can use the above equation to ob- 
tain the behavior of the scale factor at that point. When 
qi = —1 we have de- Sitter solutions (0 = constant) or 
static solutions (0 = 0). For q,; = we have a Milne evo- 
lution and when — 1 < qo < ; qo > we have accelerated 
and decelerated power law behaviours respectively. 

To obtain the exact solutions for the scale factor a(t) 
associated with the non -static 0^0 equilibrium points 



we need to have an expression for q in term of the com- 
pact variables. From the definition q we obtain: 

• (2D 

The non-invariant surface z; = Q\ is the transition 
surface between accelerated and decelerated expansions 
phases (see Fig 1). By substituting [5T] in [501 we obtain: 

e = -J(2-^)6 2 , (22) 

where Q ^ 0. The evolution of the scale factor can now 
be given directly by integrating equation 1221 

a(t) = ao(t - t a f\ (23) 

where 

A = (2" I?)- ( 24 ) 

The constants of integration can be obtained by substi- 
tuting the solutions into the original equations. As ex- 
plained in [2^|, these solutions must satisfy all the cos- 
mological equations in order to be considered physical. 

By looking at Table 1 we can distinguish two classes 
of non -static solutions. The first class {A±, D±, E±, F± 
and N±} contain solutions that are independent of n and 
lj. The fixed point B represents a static phase as men- 
tioned earlier, the points A, D and F are radiation- like 
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phases. The expanding version of the point A is a sad- 
dle for z > and repeller for z = and the other two 
points D an F are saddles. The fixed point N represent a 
matter phase and the expanding version of this point is 
a spiral - . 

The evolution of the scale factor a(t) for the fixed point 
L is a function of (n, uj, t) and for fixed point / is function 
of (n, t). The dependence of these solutions on n and/or 
uj provide us with additional degrees of freedom that can 
lead to interesting cosmological scenarios. 

When uj = 0, the fixed point L merges with N for 
n = 1, and it merges with the point D± for n = 3/4. 
When uj = 1/3 it merges with D± for n = 1, and for 
n = 4/3 it corresponds to the matter point N. In the 
case uj = or uj = 1/3, we find that for all values n for 
which this point is physical, the expansion (contraction) 
is never accelerating. 

As mentioned earlier, the evolution of the scale factor 
for the point / is independent of the equation of state 
parameter uj. For n = 5/4 the fixed point / merges with 
A, for n = 2 it merges with point E and for n = 7/12 + 
\/73/12 it is a matter point. 

We also find that for this fixed point the expansion 
(contraction) is accelerating for n > 1/2(1 + a/3)- The 
existence of this accelerated phase together with the fact 
that for n > 1/2(1 + V3) the point L is a matter-like 
point, leads to the possibility of an extremely interesting 
cosmological scenario, where it is possible in principle to 
find an orbit that starts close to the unstable accelerating 
phase /+ , evolves past the unstable radiation - like point 
D + , followed by the unstable matter point L + and finally 
ends up at the de - Sitter attractor E + . 

In Figs. 3 and 4 we have plotted two interesting orbits. 
The orbit in Fig 3 is for uj = and n = 5/4. It begins 
near the radiation- like points A+/I+, passes nearby the 
radiation - like point followed by the standard matter 
point L + and ends up at the de-Sitter attractor E + . The 
orbit in Fig 4 is for uj = -11/18 + 773/18 and n = 
7/12 + V73/12. It begins near the radiation - like point 
A + , passes nearby the matter point /+ , then close to the 
matter point L+ and ends up at the de Sitter attractor 
E+. 

It is interesting to see that for all fixed points whose 
x- coordinate goes to zero as n — > 0, are unstable. As the 
ACDM subspace lies on the x — surface and this sur- 
face is not an invariant sub- manifold, this implies that 
the behavior of the solutions of a fourth order theory 
which is close to ACDM, may be completely different 
from those of ACDM. Furthermore, this suggests that the 
best fit model to the current observational data within 
the complete state space of R + aR n , may be given by a 
non-infinitesimal value of n. 
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FIG. 1. plot of the invariant subspace z = for uj = ;n = 
5/4. The left half of the state space corresponds to collapsing 
models, while the right half contain expanding models. This is 
indicated by the subscripts of the various equilibrium points. 
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FIG. 2. For this orbit uj = and n = 5/4. 



CONCLUSION 

In this work, we have presented a careful analysis of 
the state space of the class of R + aR n theorys of gravity, 
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FIG. 3. For this orbit u = -11/18 + y/73/18 and n = 7/12 + 
V73/12. 



focusing on the R > sector with K = 0, together with 
the no-ghost condition f(R); f'(R) > 0. 

Due to the complexity of this class of gravity theories, 
the standard Hubble normalization does not lead to a 
compact dynamical variables. In order to construct vari- 
ables defining a compact dynamical system one has to use 
an appropriate normalization. In this paper we used the 
same formalism used in {28], where we absorbed all the 
negative contributions of the Friedmann equation into 
the normalization. First of all we obtained the follow- 
ing important result: For all well defined function f{R), 
with /' > and f'/Rf" invertible in terms of the dy- 
namical variables defined by (14), the FLRW universes 
with non-negative Ricci Scalar, continues to be so both 
in the future and in the past. Also an R = universe 
can never undergo a bounce in the future or past. 

Our compact analysis shows that there are more equi- 
librium points than in the corresponding non-compact 
analysis in 26]. In particular we find a new finite fixed 
point N±. Because of it's very special location in the 
phase space, it is quite difficult to obtain this point using 
the standard techniques. This point is found to represent 
a matter phase and the expanding version of this point 
is spiral^. 

Furthermore, we find that for n > 1/2(1 + ^/3) the 
phase space of R + aR n , contains two accelerated fixed 
points E+ ; 1+ , together with two other saddle points (one 
represent a radiation phase D + and the other represent 
a matter-like phase £+). Although we have obtained all 



the desired fixed points and desired stability, this does 
not necessarily imply that there is an orbit connecting 
them. Due to the fact that for n > 1/2(1 + \/3), the 
two accelerated points and the matter-like point are quite 
close to each other in the phase space, which makes it dif- 
ficult to prove the existence of an orbit connecting these 
points together with the radiation - like point. But the 
presence of all these phases in the state space of R + aR n 
makes a more detailed investigation worth pursuing. 
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